<html>
<head>
<meta charset="utf-8">
<script src="../../lib/opencv.js"></script>
<script type="text/javascript" src="https://cdnjs.cloudflare.com/ajax/libs/mathjs/5.4.0/math.min.js"></script>
</head>
<body>
<h2>Harris角点</h2>
<table border="1">
<tr>
    <th>参数</th>
    <th>值</th>
    <th>备注</th>
</tr>
<tr>
    <td>sigma</td>
    <td><input type="number" value="12" id="input_sigma"/></td>
    <td>高斯核参数，数值越大则函数越平</td>
</tr>
<tr>
    <td>高斯核尺寸</td>
    <td><input type="number" value="9" id="input_gauseKernelSize"/></td>
    <td>必须是奇数</td>
</tr>
</table>
<br>

<canvas id="canvas" style="background:DarkTurquoise "> </canvas>
<canvas id="canvas2" style="background:DeepSkyBlue "> </canvas><br>
<input id="chooseImage" type="file" accept="image/*"></input> <br>
<button id="button_run" onclick="harrisCorners(canvas)">运行</button>
<br><br>

<h4>参考链接</h4>
<li><a href="https://blog.csdn.net/newthinker_wei/article/details/45603583">Harris角点检测原理分析</a></li>


<script>
/*
TODO:
1. 考虑加入一些预处理步骤；
2. 加入非极大值抑制。
*/
const canvas = document.getElementById('canvas')
const canvas2 = document.getElementById('canvas2')
const input = document.getElementById('chooseImage')

// 图像的加载
function onImageLoaded(e) {
    const ctx = canvas.getContext('2d')
    console.log('图像加载完成')
    const image = e.target
    canvas.width = image.width
    canvas.height = image.height
    ctx.drawImage(image, 0, 0)
}

input.addEventListener('change', function(e) {
    const fileReader = new FileReader()
    fileReader.onload = function(e) {
        console.log('文件读取完成')
        const imgFile = e.target.result
        const image = new Image()
        image.onload = onImageLoaded
        image.src = imgFile
    }    
    fileReader.readAsDataURL(input.files[0])
})


function assert(value, msg) {
    if(!value) {
        if (msg)
            alert('程序出错：' + msg)
        else
            alert('程序遇到问题，停止运行！')
        throw new Object('assertion failed!')
    }
}


// 求Harris角点
function harrisCorners() {
    // 读取参数并检查参数是否合法
    const sigma = parseInt(document.getElementById('input_sigma').value)
    assert(sigma > 0, "sigma <= 0 !")
    const kernelSize = parseInt(document.getElementById('input_gauseKernelSize').value)
    assert(kernelSize % 2 == 1, "kernel Size 必须是奇数")

    // 转为灰度图
    const img = cv.imread(canvas)
    const imgGray = new cv.Mat()
    assert(img.channels() === 4)
    cv.cvtColor(img, imgGray, cv.COLOR_BGRA2GRAY)
    const rows = imgGray.rows
    const cols = imgGray.cols
    cv.imshow(canvas2, imgGray)

    // 计算x和y方向上的导数
    const dx = new cv.Mat()
    const dy = new cv.Mat()
    cv.Sobel(imgGray, dx, cv.CV_32FC1, 1, 0, 1, 1, 0, cv.BORDER_CONSTANT | 0)
    cv.Sobel(imgGray, dy, cv.CV_32FC1, 0, 1, 1, 1, 0, cv.BORDER_CONSTANT | 0)
    cv.imshow(canvas2, dy)

    // 计算一些中间变量
    const dx2 = new cv.Mat()
    cv.multiply(dx, dx, dx2)
    const dy2 = new cv.Mat()
    cv.multiply(dy, dy, dy2)
    const dxdy = new cv.Mat()
    cv.multiply(dx, dy, dxdy)
    dx2.convertTo(dx2, cv.CV_32FC1, 1 / 65535, 0)
    dy2.convertTo(dy2, cv.CV_32FC1, 1 / 65535, 0)
    dxdy.convertTo(dxdy, cv.CV_32FC1, 1 / 65535, 0)
    console.log(cv.minMaxLoc(dx))
    console.log(cv.minMaxLoc(dy))
    console.log(cv.minMaxLoc(dxdy))
    console.log('minmax dx2', cv.minMaxLoc(dx2))
    console.log('minmax dy2', cv.minMaxLoc(dy2))
    console.log('minmax dxdy', cv.minMaxLoc(dxdy))
    const matA = new cv.Mat(), matB = new cv.Mat(), matC = new cv.Mat()
    const kernel = gauseKernel(kernelSize, 1, sigma)
    cv.filter2D(dx2, matA, cv.CV_32F, kernel, new cv.Point(-1, -1), 0, cv.BORDER_CONSTANT | 0)
    cv.filter2D(dy2, matB, cv.CV_32F, kernel, new cv.Point(-1, -1), 0, cv.BORDER_CONSTANT | 0)
    cv.filter2D(dxdy, matC, cv.CV_32F, kernel, new cv.Point(-1, -1), 0, cv.BORDER_CONSTANT | 0)
    console.log('minmax A', cv.minMaxLoc(matA))
    console.log('minmax B', cv.minMaxLoc(matB))
    console.log('minmax C', cv.minMaxLoc(matC))
    imshow(canvas2, matB, {'maximumContrast': false})


    // 
    const output = cv.Mat.zeros(rows, cols, cv.CV_32FC3)
    const k = 0.1;
    const thresh = 1e-6;
    detlist = []
    trlist = []
    deltalist = []
    for (let row = 0; row < rows; row++) {
        for (let col = 0; col < cols; col++) {
            const a = matA.floatPtr(row, col)[0]
            const b = matB.floatPtr(row, col)[0]
            const c = matC.floatPtr(row, col)[0]
            det = a * b - c * c
            tr = a + b
            // r = det - k * tr * tr
            r = a * b - c * c - k * (a * a + b * b + 2 * a * b);
            detlist.push(det)
            trlist.push(tr)
            deltalist.push(r)
            pixel = output.floatPtr(row, col)
            if(r > thresh) {
                // 检测到角点
                pixel[0] = 255
            }
            else if (r < -thresh) {
                // 检测到边缘点
                pixel[1] = 255
            } 
        }
    }
    console.log("mean det: %f", math.mean(detlist))
    console.log("mean tr: %f", math.mean(tr))
    console.log("mean r: %f", math.mean(r))
    cv.imshow(canvas2, output)
}


// 计算高斯核。
// size 核尺寸，必须是单数
// a 幅度
// sigma 该值越大，核函数越是平滑
function gauseKernel(size, a, sigma) {
    assert(size % 2 == 1)
    assert(size > 0)
    const center = ~~(size / 2)
    const sigma2 = sigma**2
    const kernel = new cv.Mat(size, size, cv.CV_32F)
    for (let row = 0; row < size; row++) {
        for (let col = 0; col < size; col++) {
            kernel.floatPtr(row, col)[0] = a * Math.exp(- ((row - center)**2 / sigma2) - ((col - center)**2 / sigma2))
        }
    }
    cv.normalize(kernel, kernel, 1, 0, cv.NORM_L1, cv.CV_32F, new cv.Mat())
    let sum = 0;
    for (let row = 0; row < size; row++) {
        for (let col = 0; col < size; col++) {
            sum += kernel.floatPtr(row, col)[0]
        }
    }
    assert(Math.abs(sum - 1) < 1e-5)
    return kernel
}


function imshow(c, image, options) {
    assert(image.channels() == 1)
    const image32f = new cv.Mat()
    const minMaxLoc = cv.minMaxLoc(image);
    const min = minMaxLoc['minVal']
    const max = minMaxLoc['maxVal']
    if(options['maximumContrast'] === true && max > min) {
        const scale = 1 / (max - min) * 255
        cv.convertScaleAbs(image, image32f, scale, - min * scale)
    }
    else {
        image.convertTo(image32f, cv.CV_32F, 1, 0)
    }
    console.log(cv.minMaxLoc(image32f))
    cv.imshow(c, image32f)
}
</script>
</body>